######################### SENSITIVITY ANALYSIS ############################
####################################################################################

#Name of code file: sensitivity_analysis.R
#Purpose: Sensitivity Analysis
#Data In: network_tolerance_processed_data.csv
#Data Out: Figure C1; Tables: C1 & C2
####################################################################################

#Load packages
library(readr)
library(sensemakr)
library(xtable)

#Set Working Directory
setwd("~/Dropbox/egypt_tolerance_wp_replication/")

#Read in Data
time_data<-read_csv("data/dynamic_network_tolerance_processed_data.csv")


##############################
#FIGURE C1 and TABLES C1 & C2#
##############################

intol_change_model <- lm(intolerant_may16 ~ elite_diversity_may15 +intolerant_may15_6m + nonelite_diversity_may15 +log_total_friends +log_elite +log_time + islamist +rel_count_change_6m, data=time_data)

sensitivity <- sensemakr(model = intol_change_model, 
                         treatment = "elite_diversity_may15",
                         benchmark_covariates = c("islamist", "log_total_friends", "log_time", "nonelite_diversity_may15"),
                         kd = 1:3)

#Table C1
elite_sens<-as.data.frame(sensitivity$bounds)
xtable(elite_sens, digits=4)

pdf(file="plots/figureC1_elite.pdf", width = 11, height = 7)
plot(sensitivity, label.text=FALSE, lim=.11, lim.y=.11, cex.lab=1.5, cex.axis=1.5)
dev.off()

sensitivity <- sensemakr(model = intol_change_model, 
                         treatment = "nonelite_diversity_may15",
                         benchmark_covariates = c("islamist", "log_total_friends", "log_time", "elite_diversity_may15"),
                         kd = 1:3)


#Same with  Moderates
intol_change_model <- lm(intolerant_may16 ~ elite_diversity_may15 +intolerant_may15_6m + nonelite_diversity_mod_may15 +log_total_friends +log_elite +log_time + islamist +rel_count_change_6m, data=time_data)

sensitivity <- sensemakr(model = intol_change_model, 
                         treatment = "nonelite_diversity_mod_may15",
                         benchmark_covariates = c("islamist", "log_total_friends", "log_time", "elite_diversity_may15"),
                         kd = 1:3)

pdf(file="plots/figureC1_moderates.pdf", width = 11, height = 7)
plot(sensitivity, label.text=FALSE, cex.lab=1.5, cex.axis=1.5)
dev.off()

#Table C2
nonelite_sens_mod<-as.data.frame(sensitivity$bounds)
xtable(nonelite_sens_mod, digits=4)
